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ABSTRACT 

We  present  the  first  multi-wavelength,  high-contrast  imaging  study  confirming  the  protoplanet  em¬ 
bedded  in  the  disk  around  the  Herbig  Ae/Be  star  HD100546.  The  object  is  detected  at  L'  (~  3.8  fim) 
and  M'  (~  4.8  Rm),  but  not  at  Ks  (~  2.1  //to),  and  the  emission  consists  of  a  point  source  component 
surrounded  by  spatially  resolved  emission.  For  the  point  source  component  we  derive  apparent  mag¬ 
nitudes  of  L'  =  13.92±0.10  mag,  M'  =  13.33±0.16  mag,  and  Ks  >  15.43±0.11  mag  (3er  limit),  and  a 
separation  and  position  angle  of  (0.457±0.014)"  and  (8.4±1.4)°,  and  (0.472 ±0.014)"  and  (9.2±1.4)° 
in  L'  and  M' ,  respectively.  We  demonstrate  that  the  object  is  co-moving  with  HD100546  and  can 
reject  any  (sub-)stellar  fore-/background  object.  Fitting  a  single  temperature  blackbody  to  the  ob¬ 
served  fluxes  of  the  point  source  component  yields  an  effective  temperature  oiTeff  =  932^02  K  and  a 
radius  for  the  emitting  area  of  R  =  6.9+lg  Rjupiter-  The  best-fit  luminosity  is  L  =  (2.3+q;®)  •  10-4  Lq. 

We  quantitatively  compare  our  findings  with  predictions  from  evolutionary  and  atmospheric  models 
for  young,  gas  giant  planets,  discuss  the  possible  existence  of  a  warm,  circumplanetary  disk,  and  note 
that  the  de-projected  physical  separation  from  the  host  star  of  (53  ±  2)  au  poses  a  challenge  standard 
planet  formation  theories.  Considering  the  suspected  existence  of  an  additional  planet  orbiting  at 
~13-14  au,  HD100546  appears  to  be  an  unprecedented  laboratory  to  study  the  formation  of  multiple 
gas  giant  planets  empirically. 

Subject  headings:  planets  and  satellites:  detection  -  planets  and  satellites:  formation  -  planets  and 
satellites:  gaseous  planets  protoplanetary  disks  -  planet-disk  interactions  -  stars: 
pre-main  sequence 


1.  INTRODUCTION 

The  formation  of  gas  giant  planets  within  dust  and 
gas  rich  disks  surrounding  young  stars  is  complex  and 
theoretical  models  describing  the  immediate  formation 
process  are  barely  constrained  by  empirical  data.  At 
the  moment  there  are  two  main  theories:  the  core  accre¬ 
tion  model  based  on  the  physics  of  mutual  collisions  and 
growth  of  solid  bodies  followed  by  accretion  of  a  gaseous 
envelope  (e.g.,  Pollack  et  al.  1996),  and  the  gravitational 
instability  model  predicting  direct,  local  collapse  in  the 
outer  regions  of  circumstellar  disks  (e.g.,  Boss  2001).  Up 
to  now,  these  theories  are  indirectly  constrained,  e.g.,  by 
studying  the  chemical  and  physical  conditions  in  circum¬ 
stellar  disks  to  estimate  the  global  initial  conditions  for 
planet  formation  (e.g.,  Dutrey  et  al.  2014),  or  by  study¬ 
ing  the  composition  of  extrasolar  planets,  both  from  esti¬ 
mates  of  the  bulk  density  and  from  atmospheric  charac¬ 
terization,  to  decipher  their  formation  history  and  evo¬ 
lution  (e.g.,  Marcy  et  al.  2014;  Konopacky  et  al.  2013). 
Furthermore,  models  describing  the  luminosity  evolution 
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of  young  gas  giant  planets  are  unconstrained  by  empirical 
data  at  very  early  stages.  The  initial  specific  entropy  of 
the  objects,  which  is  dictated  by  details  of  the  gas  accre¬ 
tion  process,  is  treated  as  a  free  parameter  even  though 
its  value  has  a  significant  impact  on  the  objects  luminos¬ 
ity  in  the  first  few  hundred  million  years  (Marley  et  al. 
2007;  Spiegel  &  Burrows  2012). 

Up  to  now,  young  planet  candidates  have  been  discov¬ 
ered  inside  large  gaps  in  circumstellar  disks  surrounding 
a  few  young  stars  (Kraus  &  Ireland  2012;  Reggiani  et  al. 
2014;  Biller  et  al.  2014).  Located  within  20  au  from  their 
hosts,  they  yield  first  luminosity  estimates  of  young  gas 
giant  planets  and  suggest  that  near  the  end  of  their  for¬ 
mation,  giant  planets  have  cleared  disk  gaps  predicted 
by  theory  (e.g.,  Crida  et  al.  2006).  Until  now,  no  pro¬ 
toplanet  still  embedded  in  the  circumstellar  disk  from 
which  it  is  forming  has  been  confirmed.  We  detected  a 
candidate  protoplanet  around  the  young  star  HD100546, 
but  the  single- wavelength  data  did  not  permit  characteri¬ 
zation  nor  was  it  unambiguously  shown  to  be  a  true  com¬ 
panion  (Quanz  et  al.  2013).  Recently,  Currie  et  al.  (2014) 
recovered  this  object,  but  also  this  study  was  based  on 
single  wavelength  data  and  could  not  confirm  common 
proper  motion. 

We  now  confirm  that  this  object  is  bound  to  the  cen¬ 
tral  star  and  that  its  multi-band  photometry  is  best  ex¬ 
plained  as  a  newly  forming  gas  giant  planet  embedded  in 
the  circumstellar  disk  around  the  young  star  HD100546. 
HD100546  has  a  spectral  type  of  B9Vne  (Levenhagen  & 
Leister  2006)  and  lies  at  a  distance  of  97  ±  4  pc  (van 
Leeuwen  2007).  Based  on  photometric  measurements 
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TABLE  1 

Summary  of  observations. 


Date 

Object 

Filter  DITa 

#  of 

#  of  Airmass 

Parallactic 

Detector 

Pixel  scale 

data 

frames 

angle 

window 

cubes 

per  cube 

start  /  end 

High-contrast,  pupil  stabilized  observations 


2013/4/18 

HD100546 

L' 

0.15  s 

288 

200 

1.66-1.43 

-55.83°  /  +42.36° 

512x512  px 

27  mas/px 

HD100546 

Ks 

0.5  s 

90 

100 

1.57-1.42 

-46.34°  /  +31.95° 

512x512  px 

13  mas/px 

2013/4/19 

HD100546 

M' 

0.04  s 

244 

500 

1.61-1.43 

-50.75°  /  +24.64° 

256x256  px 

27  mas/px 

HD100546 

Ks 

0.5  s 

120 

100 

1.54-1.43 

-41.39°  /  +32.07° 

512x512  px 

13  mas/px 

Photometric  calibration  observations 


2013/4/18 

HD100546 

L' 

0.025  s 

8 

1200 

1.56-1.57 

-/- 

256x256  px 

27  mas/px 

HR6572 

L' 

0.05  s 

6 

600 

1.42-1.40 

-/- 

256x256  px 

27  mas/px 

HD100546 

Ks 

0.05  s 

4 

750 

1.58-1.59 

-/- 

256x256  px 

13  mas/px 

HR6572 

Ks 

0.05  s 

4 

750 

1.38-1.37 

-/- 

256x256  px 

13  mas/px 

2013/4/19 

HD100546 

M' 

0.01  s 

8 

2000 

~1.50 

-/- 

130x130  px 

27  mas/px 

HR6572 

M' 

0.02  s 

16 

1000 

1.71-1.61 

-/- 

130x130  px 

27  mas/px 

aDetector  integration  time,  i.e.,  exposure  time. 

covering  wavelengths  from  the  UV  to  the  infrared  van 
den  Ancker  et  al.  (1997)  estimated  the  luminosity  and 
effective  temperature  of  HD100546  to  be  L  ss  32  Lq  and 
Teg  ss  10500  K  from  which,  via  comparison  with  stellar 
evolutionary  models,  they  derived  an  age  of  >10  Myr 
and  a  mass  of  2.4  ±  0.1  M0.  A  slightly  younger  age  of 
5-10  Myr  was  estimated  by  Guimaraes  et  al.  (2006) 
from  high-resolution  optical  spectra.  These  authors  used 
the  spectra  to  infer  the  effective  temperature  and  sur¬ 
face  gravity  of  HD100546  and  then  compared  the  values 
to  stellar  evolutionary  tracks  to  get  an  age.  As  the  age 
of  the  star  is  an  important  parameter  in  the  context  of 
(gas  giant)  planet  formation,  we  will  consider  a  range  of 
ages  between  1-10  Myr  throughout  this  paper,  where 
the  youngest  age  is  motivated  by  the  idea  that  our  ob¬ 
ject  might  be  younger  than  the  star  as  it  is  still  in  the 
process  of  formation.  The  large  (r  >  300  au)  gas  and 
dust  disk  around  HD  100546  has  been  spatially  resolved 
in  scattered  light  as  well  as  in  thermal  emission  (Pantin 
et  al.  2000;  Augereau  et  al.  2001;  Grady  et  al.  2001;  Liu 
et  al.  2003;  Leinert  et  al.  2004;  Ardila  et  al.  2007;  Quanz 
et  al.  2011;  Avenhaus  et  al.  2014;  Pineda  et  al.  2014; 
Walsh  et  al.  2014).  The  star  is  actively  accreting  mate¬ 
rial  from  the  innermost  disk  regions  (e.g.,  Deleuil  et  al. 
2004;  Grady  et  al.  2005;  Guimaraes  et  al.  2006),  but  there 
is  observational  evidence  for  a  disk  gap  stretching  from 
a  few  au  out  to  roughly  14  au  that  is  possibly  created 
by  an  orbiting  companion  (Bouwman  et  al.  2003;  Grady 
et  al.  2005;  Acke  &  van  den  Ancker  2006;  Benisty  et  al. 
2010;  Quanz  et  al.  2011;  Avenhaus  et  al.  2014;  Brittain 
et  al.  2013,  2014).  This  inner  companion  together  with 
the  newly  forming  companion  in  the  outer  disk  discussed 
in  the  following  makes  HD100546  currently  the  best  can¬ 
didate  for  a  system,  where  we  can  directly  study  the 
formation  of  multiple  planets  and  their  interaction  with 
the  circumstellar  disk. 

2.  OBSERVATIONS 

We  re-observed  HD100546  on  April  18  and  19,  2013, 
with  the  NACO  instrument  (Lenzen  et  al.  2003;  Rousset 
et  al.  2003)  installed  at  one  of  the  8.2-m  Utility  Tele¬ 
scopes  of  the  Very  Large  Telescope  at  the  Paranal  Ob¬ 
servatory  of  the  European  Southern  Observatory.  Data 
were  taken  in  the  L'  (Aes  =  3.770  /im)  and  Ks  (Aeff  = 


2.124  fim)  filters  in  night  1,  and  in  the  M'  (Aeff  = 
4.755  pirn)  and  again  Ks  filters  in  night  2.  The  obser¬ 
vations  were  carried  out  in  pupil  tracking  mode  leading 
to  a  natural  rotation  of  the  cameras  field  of  view  fol¬ 
lowing  the  changes  in  parallactic  angle  over  the  course 
of  the  observing  sequence.  Each  night  we  switched  be¬ 
tween  the  different  filters  several  times  to  ensure  com¬ 
parable  field  rotation  in  all  datasets.  To  correct  for  bad 
detector  pixels  and  allow  for  subtraction  of  thermal  back¬ 
ground  emission,  the  objects  were  moved  to  different  po¬ 
sitions  on  the  detector  every  30  to  60  seconds  in  the  L' 
and  M'  filter  and  roughly  every  90  seconds  in  the  Ks 
filter.  For  the  data  taken  in  the  L'  filter  the  Apodiz- 
ing  Phase  Plate  (APP)  coronagraph  (Kenworthy  et  al. 
2010;  Quanz  et  al.  2010)  was  used  to  enhance  the  con¬ 
trast  performance  of  the  instrument.  No  coronagraph 
was  used  for  the  M'  and  Ks  filter  observations.  In  the 
high-contrast  datasets,  in  order  to  increase  the  sensitivity 
to  faint  companions,  the  central  few  pixels  of  the  stellar 
Point  Spread  Function  (PSF)  were  saturated,  resulting 
in  the  need  for  additional,  unsaturated  datasets  to  de¬ 
termine  the  exact  photometry  of  HD100546  and  thus  its 
companion.  As  the  observing  conditions  were  photomet¬ 
ric  in  both  nights,  the  photometric  standard  star  HR6572 
was  observed  as  reference  target.  The  observations  are 
summarized  in  Table  1. 

3.  DATA  REDUCTION 
3.1.  Basic  steps 

Pupil-stabilized,  data  Ks  filter:  Individual  frames  were 
flat-fielded,  dark-subtracted  and  bad  pixel  corrected  (5-cr 
clipping  in  9  pixel  box).  For  image  alignment,  a  Moffat 
profile  was  fitted  to  a  reference  image  and  to  the  individ¬ 
ual  science  images  to  determine  the  spatial  offset  between 
them.  The  individual  images  were  then  scaled  up  by  a 
factor  of  10,  shifted  by  the  offset  and  then  scaled  back  to 
the  original  image  size. 

Pupil-stabilized  data  L'  ’filter  with  APP:  Individual  ex¬ 
posures  from  subsequent  detector  positions  were  pairwise 
subtracted  (to  subtract  background  emission  and  dark 
current),  and  bad  pixel  corrected  (5-cr  clipping  in  9  pixel 
box).  For  image  alignment,  given  the  strong  asymmetry 
of  the  APP  PSF  with  the  diffraction  rings  being  sup¬ 
pressed  on  one  side  of  the  central  star  (Codona  et  al. 
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2006),  we  cross-correlated  the  individual  science  images 
with  a  reference  image  to  determine  the  spatial  offset  be¬ 
tween  them.  The  individual  images  were  then  scaled  up 
by  a  factor  of  10,  shifted  by  the  offset  and  then  scaled 
back  to  twice  the  original  image  size. 

Pupil-stabilized  M'  filter:  Individual  frames  were  dark- 
subtracted  and  bad  pixels  were  flagged  with  a  bad  pixel 
mask.  Background  emission  was  then  subtracted  using  a 
novel  principal  component-based  approach  (Quanz  et  al. 
in  prep):  The  star  followed  a  4-point  dither  pattern  on 
the  detector,  where  each  position  was  roughly  centered 
in  each  of  the  4  detector  quadrants.  Hence,  throughout 
the  full  stack  of  image  cubes,  for  a  given  quadrant,  the 
star  was  not  located  in  this  quadrant  ~75%  of  the  time. 
These  starless  frames  were  used  for  the  decomposition 
of  the  background  emission  into  principal  components 
(PCs).  In  order  to  fit  the  background  in  a  given  quadrant 
in  those  frames  where  the  star  was  present,  the  innermost 
parts  around  the  star  were  masked  out  and  then  the  PCs 
were  fitted  to  the  remaining  area.  The  background  for 
a  given  frame  was  then  constructed  from  this  fit,  which 
automatically  interpolated  those  innermost  regions  con¬ 
taining  the  star  that  had  been  masked  out  before.  For 
image  alignment,  a  Moffat  profile  was  fitted  to  a  reference 
image  and  to  the  individual  science  images  to  compute 
the  offset  between  them.  The  individual  images  were 
then  scaled  up  by  a  factor  of  10,  shifted  by  the  offset  and 
then  scaled  back  to  twice  the  original  image  size. 

For  all  filters,  individual  images  with  poor  adaptive- 
optics  correction  were  disregarded  from  any  subsequent 
analyses.  Those  images  were  identified  during  the  initial 
data  reduction  steps  by  comparing  the  PSF  peak  flux  of 
consecutive  images:  in  case  the  flux  dropped  by  ~40% 
from  one  image  to  the  next  we  flagged  and  visually  in¬ 
spected  the  image.  In  total,  less  than  1%  of  the  images 
suffered  from  bad  AO  performance  and  were  rejected. 

3.2.  PSF  subtraction 

To  reveal  the  existence  of  possible  faint  nearby  com¬ 
panions,  the  stellar  PSF  is  subtracted  from  the  individual 
images.  This  was  achieved  using  the  PynPoint  data 
analysis  package  (Amara  &  Quanz  2012;  Arnara  et  al. 
2014)  that  is  based  on  a  PC  analysis  algorithm  to  model 
the  light  distribution  of  the  PSF.  The  PCs  are  fitted  to 
each  individual  image  and  the  stellar  PSF  is  subtracted. 
All  images  are  then  rotated  to  a  common  sky  orientation, 
averaged  and  convolved  with  a  Gaussian  kernel  with  a 
full-width-half-maximum  of  half  the  size  of  the  stellar 
PSF  to  increase  the  signal-to-noise  (SNR)  of  any  faint 
companion.  The  size  of  the  stellar  PSF  was  determined 
from  the  average  of  the  unsaturated  data  sets  in  each 
filter. 

4.  OBSERVATIONAL  RESULTS 
4.1.  Recovery  of  the  protoplanet  in  L'  and  M' 

The  fully  processed  L'  and  M'  images  reveal  an  emis¬ 
sion  source  North  of  the  star.  The  source  is  spatially  re¬ 
solved  in  both  images,  but  one  dimensional  cuts  through 
the  object’s  flux  distribution  in  the  radial  direction  re¬ 
veal  a  point  source  component  that  dominates  the  flux. 
Hence,  the  object  is  best  explained  by  a  spatially  un¬ 
resolved  point  source  component  surrounded  by  a  spa¬ 
tially  resolved  extended  emission  component  (Figure  1, 


top  row) .  The  total  SNR  of  the  detection  was  estimated 
following  the  prescription  of  Mawet  et  al.  (2014)  and 
their  equation  9,  taking  into  account  the  peak  flux  of 
the  planet  and  the  standard  deviation  in  a  number  of 
background  pixels  each  representing  a  statistically  inde¬ 
pendent  resolution  element  at  the  separation  from  the 
star  that  is  of  interest.  As  the  final  PynPoint  images 
were  convolved  with  a  Gaussian  filter,  individual  pixels 
separated  enough  to  be  statistically  independent,  i.e. ,  one 
pixel  per  resolution  element,  were  chosen  as  background 
pixels.  The  background  pixels  were  selected  from  2  con¬ 
centric  rings  around  the  star  with  radii  of  ~0.46"  and 
~0.54",  respectively,  to  have  a  comparable  separation 
from  the  central  star  as  the  object  and  the  surrounding 
extended  emission  component.  For  the  L'  images  taken 
with  the  APP  coronagraph  only  the  high-contrast  side  of 
the  APP  was  considered  for  selecting  background  pixels; 
for  the  M'  images  complete  circles  were  considered.  Ex¬ 
cluded  were  the  immediate  surroundings  of  the  source  as 
well  as  the  extended  emission  region  east  of  the  star  (see, 
bottom  row  in  Figure  1  and  section  6).  This  resulted  in 
19  and  32  background  pixels  for  the  L'  and  M'  images, 
respectively.  We  then  computed  the  SNR  for  a  grid  of 
different  reduction  parameters  varying: 

1.  The  number  of  individual  images  that  were  aver¬ 
aged  before  PynPoint  was  run  (5,  10  and  50  in 
L'\  10,  20  and  100  in  M') 

2.  The  number  of  PCs  used  to  fit  the  PSF  (between 
10  and  50  in  M'  and  between  10  and  120  in  L', 
always  in  steps  of  10) 

3.  The  exact  location  of  the  reference  pixels  (shifting 
all  pixels  simultaneously  by  ±1  pixel  in  x  and  y 
resulting  in  a  total  number  of  9  sets  of  reference 
pixels) 

In  all  of  these  reductions  the  SNR  of  the  object  was 
always  >4  in  L'  and  >3  in  M'  and  the  mean  SNR  in 
L'  and  M'  was  «7.9  and  «4.6.  Irrespective  of  the  exact 
location  of  the  reference  pixels,  the  highest  SNR  values 
were  always  achieved  for  70  PCs  and  pre-averaging  of 
10  images  in  L'  filter  (average  SNRssll.4)6  and  for  10 
PCs  and  pre-averaging  over  100  images  in  the  M'  filter 
(average  SNR«6.6).  The  resulting  images  based  on  these 
parameter  were  used  for  all  subsequent  analyses  and  are 
shown  in  Figure  1. 

We  note  that  given  the  small  number  of  background 
pixels  it  is  not  possible  to  robustly  constrain  their  under¬ 
lying  distribution  without  making  some  ad-hoc  assump¬ 
tion  of  their  general  functional  form.  In  consequence, 
without  assuming  a  functional  form  for  the  underlying 
distribution,  we  cannot  assign  a  confidence  level  to  our 
detection.  However,  given  that  the  object  was  already 
detected  previously  in  independent  datasets  (Quanz  et  al. 
2013;  Currie  et  al.  2014)  and  is  now  re-detected  in  two 
different  filters  makes  us  believe  that  a  statistical  outlier 
or  instrumental  artifact  is  extremely  unlikely. 

6  We  note  that  the  SNR  of  the  L'  detection  in  the  discovery 
paper  (SNRR315;  Quanz  et  al.  2013)  was  overstimated  as  we  did  not 
carry  out  a  rigorous  SNR  assessment  as  we  did  here  (incl.  varying 
in  the  reduction  parameters  and  background  pixels,  correcting  for 
any  offset  in  the  mean  level  of  the  background  pixels,  etc.) 
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Fig.  1. —  Top  row:  Final  PSF-subtracted  images  of  the  vicinity  around  HD100546  in  the  L'  filter  (left  panel)  and  M'  filter  (right  panel). 
The  dark  spot  in  the  image  center  indicates  the  location  of  the  central  star.  The  innermost  regions  have  been  masked  out  during  the  data 
analysis  as  they  are  dominated  by  subtraction  residuals.  Bottom  row:  Same  as  above  but  after  subtraction  of  the  point  source  component. 
Residual  resolved  structures  remain  visible  in  the  vicinity  of  the  protoplanet.  Additional  extended  emission  is  present  to  the  southeast  of 
the  star  mainly  in  the  M'  filter  (see  also  Section  6).  North  is  up  and  East  to  the  left  in  all  images. 


4.2.  Contrast  and  astrometry  of  the  protoplanet  in  L' 

and  M' 

The  PSF  subtraction  step  affects  the  exact  location 
and  brightness  of  any  companion.  To  estimate  the  pro¬ 
toplanet’s  contrast  and  position,  artificial  negative  ob¬ 
jects  (covering  a  grid  of  varying  brightness  (in  steps  of 
0.1  mag)  and  location  (in  steps  of  0.25  pixels  in  x  and  y, 
respectively))  were  inserted  in  the  individual  input  im¬ 
ages  and  PynPoint  was  re-run.  An  unsaturated  PSF- 
core  from  the  photometric  datasets  was  used  as  template 


to  create  the  fake  sources.  To  determine  their  flux  lev¬ 
els  in  the  saturated  datasets  differences  in  airmass  dur¬ 
ing  the  observations  and  the  difference  in  exposure  time 
between  the  saturated  and  unsaturated  datasets  were 
taken  into  account.  To  subtract  the  stellar  PSF  we  used 
the  same  number  of  PCs  as  in  the  final  images  defined 
above.  However,  a  new  set  of  PCs  was  constructed  for 
each  new  dataset  containing  fake  native  planets.  Based 
on  the  previous  results  from  the  initial  discovery  paper 
(Quanz  et  al.  2013)  and  from  the  final  images  shown  in 
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Figure  1,  the  baseline  assumption  was  that  the  emission 
from  the  companion  consists  of  a  point  source  component 
surrounded  by  spatially  resolved,  extended  emission  com¬ 
ponent.  To  determine  the  contrast  and  location  of  the 
point  source  component  the  goal  was  to  subtract  only 
so  much  of  the  emission  so  that  the  remaining  flux  is 
comparable  to  the  flux  level  in  the  immediate  vicinity  of 
the  object  (i.e.,  the  extended  emission  component)  with¬ 
out  creating  strong  structural  asymmetries  in  the  final 
image.  To  evaluate  the  impact  of  removing  flux  from 
the  point  source  we  quantify  the  level  of  curvature  at 
the  objects  position  using  the  determinant  of  the  Hes¬ 
sian  matrix  of  the  image  at  that  point.  In  an  analogous 
way  to  second  derivatives  in  ID,  excess  flux  from  a  com¬ 
pact  source  would  lead  to  a  positive  determinant  while  an 
over  subtraction  would  lead  to  a  hole  and  thus  a  negative 
determinant.  This  allowed  us  to  determine  the  best-fit 
subtraction,  which  was  confirmed  through  visual  inspec¬ 
tion:  the  residual  images,  once  the  flux  from  the  point 
source  had  been  removed  in  this  way,  looked  consistent 
with  a  smooth  broader  background  emission  region  (see, 
Figure  1  bottom  row).  For  further  discussion  of  image 
feature  detection  using  Hessian  matrix  methods  and  oth¬ 
ers  we  refer  the  reader  to,  e.g.,  Lindeberg  (1998)  and  Bay 
et  al.  (2008). 

For  the  L'  images  this  resulted  in  a  contrast  of  AL'  = 
9.4  ±0.1  mag  between  HD100546  and  the  compact  emis¬ 
sion  component.  The  accuracy  of  the  location  of  the  com¬ 
pact  emission  component  was  ±0.5  pixels  in  the  Pyn- 
Point  images,  which  translates  into  ±0.25  pixels  in  the 
original  image  size  (see  Section  3.1  above).  This  error 
does  not  include  systematic  uncertainties  if  we  had  cho¬ 
sen  a  different  final  image  (i.e.,  with  a  different  number 
of  PC  used  to  subtract  the  stellar  PSF,  see  above).  Fake 
sources  that  are  brighter/fainter  or  more  offset  than  the 
values  quoted  above  left  a  measurable  level  of  curvature 
and  clearly  visible  brightness  asymmetries  in  the  final 
subtracted  images.  The  same  analysis  carried  out  for  the 
M'  data  yielded  a  contrast  of  AM  =  9.2  ±  0.15  mag  be¬ 
tween  HD100546  and  the  compact  emission  component 
and  the  same  accuracy  in  the  object’s  location. 

It  is  worth  noting  that  after  subtraction  of  the  best- 
fit  negative  point  source  the  remaining  extended  emis¬ 
sion  is  not  only  elongated  in  the  radial  direction  in  the 
final  images,  but  now  also  left  and  right  of  the  point 
source  component  persisting  emission  structures  appear 
(see,  Figure  1  bottom  row).  This  can  be  understood 
if  one  keeps  in  mind  that  normally  any  reduction  algo¬ 
rithm  for  pupil  stabilized  images  creates  negative  holes 
left  and  right  of  a  point  source  as  it  tries  to  cancel  it 
out.  In  the  final  PynPoint  images  no  negative  holes 
are  present,  but  the  source  has  an  unusually  elongated 
shape,  which  we  interpret  as  PynPoint  canceling  out 
flux  left  and  right  of  the  point  source  component.  Part 
of  this  flux  reappears  after  the  point  source  component 
is  subtracted.  Without  extensive  modeling  the  precise 
shape  and  brightness  distribution  of  the  extended  emis¬ 
sion  component  cannot  be  determined,  however,  and  we 
basically  assume  that  it  is  flat  for  the  analysis  described 
here. 

The  final  astrometry  of  the  object  was  derived  using 
the  location  of  the  best-fit  negative  fake  planet  and  the 
location  of  the  central  star  (with  an  uncertainty  of  ±0.2 
pixels  in  x  and  y)  and  we  found  (0.457  ±  0.014)"  and 


(0.472  ±  0.014)"  for  the  separation  and  (8.4  ±  1.4)°  and 
(9.2  ±  1.4)°  for  the  position  angle  in  the  L'  and  M'  fil¬ 
ter,  respectively.  These  values  are  consistent  with  those 
found  in  our  discovery  paper  (Quanz  et  al.  2013).  The 
errors  are  calculated  from  error  propagation  in  r  (the 
separation  derived  from  Ax  and  Ay  between  planet  and 
star)  and  in  Ax/ Ay  (which  is  used  to  compute  the  po¬ 
sition  angle  via  PA  =  tan-1  {Ax/ Ay)).  A  systematic 
error  in  the  position  angle  from  an  unknown  true  north 
orientation  of  the  camera  field-of-view  is  not  included, 
but  is  estimated  to  be  <  1°.  Assuming  an  inclination 
and  position  angle  of  the  circumstellar  disk  of  42°  and 
145°,  respectively  (Pineda  et  al.  2014),  the  average  de- 
projected  physical  separation  of  the  compact  object  is 
(53  ±  2)  au.  The  error  assumes  the  same  uncertainties 
in  r  and  PA  as  for  the  individual  measurements  of  these 
parameters  in  L'  and  M'  band  and  denotes  the  results 
for  the  most  extreme  combinations  of  r  and  PA  allowed 
in  the  given  range.  Uncertainties  in  the  disk  inclination 
and  disk  position  angle  are  not  included. 

4.3.  Non-detection  of  the  protoplanet  in  Ks 

While  the  protoplanet  is  clearly  detected  in  the  L'  and 
M'  filter,  the  object  is  not  detected  in  any  of  the  two 
Ks  datasets  that  were  obtained.  To  estimate  an  upper 
limit  of  the  objects  I\s  brightness  we  inserted  fake  plan¬ 
ets  with  varying  brightness  (in  steps  of  0.1  mag)  at  the 
expected  location  derived  from  the  L'  and  M'  data  and 
re-ran  PynPoint.  For  this  we  used  the  Ks  dataset  from 
night  1  as  overall  the  observing  conditions  and  the  AO 
performance  were  better  compared  to  night  2.  To  com¬ 
pensate  for  a  lower  Strehl  ratio  and  higher  variability  in 
the  AO  correction  in  the  Ks  filter  compared  to  L'  or  M' , 
we  selected  the  best  50%  of  the  image  frames  from  night 
1.  For  this  we  subtracted  the  mean  Ks  image  from  night 
1  from  all  individual  images  and  computed  the  residual 
noise  in  the  resulting  image  by  determining  the  standard 
deviation  of  all  pixels.  We  then  sorted  all  images  by 
their  noise  level  and  took  the  best  50%,  i.e.,  those  with 
the  lowest  noise,  as  our  input  images  for  further  analy¬ 
ses.  To  estimate  the  upper  limit  for  the  Ks  brightness 
we  applied  the  same  approach  as  done  for  the  L'  and  M' 
detections  (see,  section  4.1;  cf.  Mawet  et  al.  2014).  We 
estimated  the  SNR  of  fake  companions  for  different  num¬ 
bers  of  PCs  and  2  different  sets  of  input  data:  one,  where 
each  of  the  selected  images  was  used,  and  one,  where  10 
consecutive  images  were  stacked  before  PynPoint  was 
run  (similar  to  the  analysis  in  L'  and  M').  As  back¬ 
ground  pixels  we  used  24  statistically  independent  pixels 
from  a  concentric  ring  around  the  central  star  with  a  ra¬ 
dius  corresponding  to  the  separation  between  star  and 
protoplanet.  We  found  that  for  a  contrast  of  AKS  =  9.6 
mag  a  fake  companion  would  have  been  detected  with 
an  average  SNR  of  «  3.7  in  both  sets  of  input  data  for 
reductions  with  5,  10,  20,  and  40  PCs  (see,  Figure  2). 
Given  that  we  know  from  the  L'  and  M'  datasets  that 
there  is  a  source  at  this  location  we  use  this  contrast  as 
our  3-<t  upper  limit. 

We  note  that  already  Boccaletti  et  al.  (2013)  did  not 
detect  the  protoplanet  in  the  Ks  filter  in  archival  data 
from  Gemini/NICI.  The  detection  limits  we  achieve  here 
are  somewhat  deeper,  further  emphasizing  the  very  red 
colors  of  the  object  (see  below). 
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Fig.  2. —  Derivation  of  the  upper  flux  limit  in  the  Ka  filter.  Left:  Final  Ka  image  reduced  with  10  principal  components  in  PynPoint. 
No  significant  point  source  is  detected.  Right:  Same  image  with  an  artificial  planet  with  A Ka  =  9.6  mag  inserted  in  the  input  data  at  the 
expected  location  of  the  companion.  The  object  is  detected  with  an  SNR  ~  3  (black  arrow).  North  is  up  and  East  to  the  left. 


4.4.  Photometry  of  HD1005f6  and  the  protoplanet 

We  observed  HR6572  (an  AOV  star)  as  photometric 
standard  star  in  both  nights  as  the  observing  conditions 
were  photometric.  Unfortunately,  HR.6572  has  no  listed 
M'  magnitude  and  hence  we  assumed  that  M'  =  L' . 
For  the  photometry,  each  data  cube  (see,  Table  1)  was 
reduced  and  analyzed  individually  including  bad  pixel 
cleaning,  background  subtraction,  image  alignment  and 
averaging  (see,  Section  3.1).  The  stellar  flux  was  mea¬ 
sured  in  the  final  image  in  an  aperture  with  3  pixels  ra¬ 
dius  in  the  Ks  and  L'  filter,  and  in  an  aperture  with  4 
pixels  radius  in  the  M'  filter.  The  average  flux  of  all  final 
images  was  taken  as  the  final  flux  and  the  standard  devi¬ 
ation  of  the  flux  of  all  final  images  divided  by  the  square 
root  of  the  number  of  final  images  was  taken  as  error 
on  the  flux  measurement  (i.e. ,  the  standard  deviation  of 
the  mean).  The  final  apparent  magnitudes  for  HD100546 
were  obtained  by  comparing  its  final  fluxes  in  the  differ¬ 
ent  filters  to  those  of  HR6572  (normalized  to  the  same 
integration  time)  and  using  the  cataloged  magnitudes  of 
HR6572  as  reference  points.  To  correct  for  the  difference 
in  airmass  between  HD100546  and  the  HR6572,  we  used 
the  Paranal  extinction  values  listed  on  the  ESO/NACO 
webpage  ( Ks :  0.07  mag,  L'\  0.08  mag,  M 0.15±0.05 
mag7  ).  The  final  photometric  errors  (in  magnitudes)  for 
HD100546  were  computed  via 

_  ( (M'  \2  ,  (  M'  \2 

°total,HD100546  ~  U(Jobs,HD100546l  '  l°obs,HR6572l 

M'  \2  ,  /  M'  \2\ 1^2 

'  \^cat,HR6572/  '  v  ciirmass/  I 

with 

7  For  the  M'  filter  ESO  provides  a  range  for  the  extinction  value 
of  0.1— 0.2  mag.  We  chose  to  use  the  midpoint  of  this  range  and 
included  an  error  term  in  subsequent  analysis  to  reflect  the  uncer¬ 
tainty. 


^obs  hdioo546:  observed  standard  deviation  of  mean  flux 
’of  HD100546 

crobsHR6572:  observed  standard  deviation  of  mean  flux 
of  HR6572 

CTcat  HR6572*  photometric  error  of  HR6572  in  catalog 
oA7  „  •  error  in  airmass  correction. 

dll  llldoo 

For  L '  and  Ks  the  magnitude  errors  were  computed 
accordingly,  with  the  only  difference  being  that  no  error 
in  the  airmass  correction  was  included.  Using  the  con¬ 
trast  measurements  from  the  previous  section  and  the 
apparent  magnitudes  for  HD100546,  we  derived  the  ap¬ 
parent  magnitudes  of  the  compact  emission  component, 
i.e.,  the  protoplanet.  For  this,  we  also  investigated  how 
color  correction  terms  due  to  the  red  color  of  the  object 
and  the  different  filter  systems  in  NACO  and  in  the  ref¬ 
erence  star  catalog  might  affect  the  results.  We  found 
this  effect  to  be  <2%,  which  is  smaller  than  the  final 
uncertainties  in  the  apparent  magnitudes  of  the  object. 
We  summarize  the  photometric  results  in  Table  2.  The 
errors  for  HR6572  comprise  the  second  and  third  term 
from  the  right  hand  side  of  the  equation  above. 

5.  ANALYSIS 

Throughout  the  paper  we  already  referred  to  the  de¬ 
tected  point  source  component  as  ’protoplanet’.  In  the 
following  we  show  that  its  observed  properties  are  indeed 
best  explained  with  a  young,  forming  planet. 

5.1.  Rejection  of  fore- /background  objects 

To  check  for  common  proper  motion  between  the  pro¬ 
toplanet  and  HD100546  the  L'  data  presented  here  we 
taken  with  exactly  the  same  observational  setup  as  the 
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TABLE  2 

Observed  and  derived  apparent  magnitudes. 


Object 

Ks  [mag] 

L'  [mag] 

M'  [mag] 

HR6572 

5.738±0.016 

5.726±0.009 

5.73a±0.014 

HD100546 

5.83±0.06 

4.52±0.02 

4.13±0.05 

HD100546  b 

>15.43b±0.06 

13.92±0.10 

13.33±0.16 

>13.59b,c±0.06 

12.75c±0.10 

12.33c±0.16 

aHR6572  (an  AOV  star)  has  no  listed  M'  magnitude  in  the  stan¬ 
dard  star  catalog,  hence  we  assume  M'  Ri  L'(r  Ka),  but  increase 
the  associated  error. 

b3(r  detection  limit;  the  error  reflects  the  uncertainty  in  the  pho¬ 
tometry  of  HD100546. 

cValues  after  correcting  for  dust  extinction  effects  (section  5.3.2.) 

data  presented  in  the  initial  discovery  paper  (Quanz  et  al. 
2013).  A  proper  motion  analysis  based  on  the  measured 
astrometry  of  the  protoplanet  and  the  parallactic  motion 
and  proper  motion  of  HD100546  between  the  two  epoch 
shows  that  the  object  is  inconsistent  with  a  stationary 
background  source  (see,  Figure  3).  The  analysis  assumes 
a  proper  motion  of  HD100546  in  right  ascension  and  dec¬ 
lination  of  —38.93  mas/yr  and  0.29  mas/yr,  respectively 
(van  Leeuwen  2007).  The  astrometric  accuracy  of  the 
data  presented  here  is  discussed  above.  A  re-analysis  of 
the  data  from  epoch  1  led  to  an  uncertainty  in  the  proto¬ 
planets  location  of  ~0.75  pixels  in  x  and  y  (inch  uncer¬ 
tainties  in  the  location  of  the  central  star).  We  made  sure 
that  both  datasets,  from  epoch  1  and  2,  were  aligned  to 
the  same  reference  image  to  minimize  systematic  offsets 
in  determining  the  location  of  the  protoplanet.  As  the 
data  presented  here  have  higher  SNR  than  the  epoch  1 
data,  they  provide  a  better  precision  in  the  protoplanets 
location. 

In  addition  to  this,  the  combination  of  apparent  L'  and 
M'  band  brightness  and  color  are  inconsistent  with  those 
of  any  stellar  foreground  or  background  object.  Only 
some  very  cool  brown  dwarfs  with  spectral  types  of  T6 
and  later  show  similar  properties,  but,  without  excep¬ 
tion,  these  objects  are  located  in  the  immediate  neigh¬ 
borhood  from  the  Sun  with  distances  typically  <10  pc 
and  hence  very  large  proper  motions  (Golimowski  et  al. 
2004;  Faherty  et  al.  2009).  Given  that  the  location  of 
our  object  has  not  changed  compared  to  data  from  2011 
we  can  exclude  that  it  is  a  nearby  object.  Rather,  the 
red  color  of  the  compact  component  combined  with  the 
extended  emission  component  suggests  that  the  object  is 
associated  with  the  circumstellar  disk  of  HD100546. 

5.2.  Rejection  of  scattered  light  from  the  circumstellar 

disk 

As  HD100546  is  surrounded  by  a  large,  flared  circum¬ 
stellar  disk  that  has  been  detected  in  scattered  light 
at  multiple  wavelengths,  the  detected  emission  (point 
source  +  extended  component)  could,  in  principle,  be 
starlight  reflected  from  the  disk’s  dusty  surface  layer. 
However,  using  polarized  light  as  tracer  for  scattered 
light,  this  seems  rather  unlikely.  No  local  brightness 
maximum  is  seen  in  high-contrast  polarized  light  images 
of  HD100546  at  the  location  of  the  object,  neither  in  the 
NIR  at  H  or  Ks  nor  in  L'  (Quanz  et  al.  2011;  Avenhaus 
et  al.  2014). 

5.3.  Effective  temperature  and  emitting  area 


R.A.  offset  (arcsec) 


Fig.  3. —  Proper  motion  analysis  based  on  the  observed  astrom¬ 
etry  of  the  protoplanet  in  the  two  epochs.  The  x  and  y  axes  show 
the  offset  with  respect  to  the  central  star.  The  red  crosses  de¬ 
note  the  location  of  the  protoplanet  relative  to  the  star  in  the  two 
epochs.  The  black  line  shows  the  expected  motion  of  a  background 
source  between  the  two  epochs  based  on  the  parallactic  motion  and 
proper  motion  of  HD100546. 

5.3.1.  Without  dust  extinction  effects 

Rejecting  scattered  light  as  origin  for  the  observed  flux, 
leaves  thermal  emission  coming  from  (within)  the  cir¬ 
cumstellar  disk  as  a  possible  source.  As  derived  above 
the  de-projected  physical  separation  of  the  compact  ob¬ 
ject  is  (53  ±  2)  au.  At  this  separation  from  the  central 
star,  radiative  transfer  models  predict  a  temperature  of 
~50  K  in  the  mid-plane  of  the  HD100546  circumstellar 
disk  (Mulders  et  al.  2011),  which  is  inconsistent  with  the 
observed  L'  —  M'  color.  A  local  extra  source  of  energy  is 
required. 

Assuming  a  distance  of  (97  ±  4)  pc,  we  used  black- 
body  emission  to  estimate  the  effective  temperature  and 
emitting  area  of  the  detected  point  source  component. 
We  computed  blackbody  fluxes  for  a  2D  grid  of  tempera¬ 
tures  Teff  (from  500  to  3000  K  in  steps  of  5  K)  and  radii 
R  (from  1  to  25  Rjupiter  hi  steps  of  0.1  Rjupiter)  and  con¬ 
volved  them  with  the  NACO  filter  transmission  curves. 
We  then  computed  a  y2-grid  (each  cell  corresponding  to 
a  certain  Teff  —  R  combination)  by  fitting  the  blackbody 
fluxes  to  the  observed  fluxes8  and  converted  this  grid 
into  a  likelihood  grid  where  the  likelihood  in  each  cell  is 
p  oc  exp(— X2/2).  In  these  fits  the  non-detection  in  the 
Ks  filter  was  explicitly  taken  into  account  by  measuring 
the  average  flux  at  the  expected  location  of  the  planet  in 
the  8  final  Ks  images  used  in  the  SNR  analysis  described 
in  section  4.3.  This  flux  was  used  as  ’observed’  flux  in 
the  Ks  filter.  In  all  the  blackbody  fits,  uncertainties  in 
the  distance  estimate  and  photometry  of  HD100546  as 
well  as  in  the  photometry  of  the  protoplanet  (in  case  of 
the  Ks  filter  the  upper  limit)  were  taken  into  account. 

The  probability  distributions  for  Te//  and  R  were  com¬ 
puted  by  marginalizing  over  the  other  parameter  in  the 
likelihood  grid  and  normalizing  the  resulting  distribu¬ 
tion.  These  were  then  used  to  compute  the  expectation 
values  and  confidence  levels  for  Teff  and  R.  This  exercise 
yielded  Teff  =  932+202  K  and  R  =  6.9+2'g  Rjupiter  for 
the  effective  temperature  and  radius  of  the  emitting  area, 

8  We  assumed  the  following  zero  points  (in  erg/cm2/s/A): 
4.501  ■  10-11,  5.151  ■  10-12  and  2.117-  10~12  for  Ks,  L' ,  and  M' , 
respectively  (see,  http://svo2.cab.inta-csic.es/theory/fps3/) 
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Fig.  4. —  Results  from  the  y2-fits  of  blackbodies  with  varying  temperature  and  size  to  our  Ks,  L'  and  M'  photometry.  2D  confidence 
levels  are  shown  in  the  bottom-left  panel  with  the  black,  green  and  red  colored  regions  denoting  specific  confidence  regions.  The  probability 
distributions  for  Teff  and  R  are  shown  in  the  top  and  bottom-right  panels,  respectively,  with  the  dashed  lines  indicating  the  mean  value 
and  the  dotted-lines  enclosing  the  ~68 %  confidence  interval. 


respectively  (the  error  bars  define  the  ~68%  confidence 
interval) . 

We  also  normalized  the  full  likelihood  grid  to  identify 
those  combinations  in  the  Teff  —  R  space  that  corre¬ 
sponded  to  certain  confidence  regions  (Figure  4).  The 
best-fit  total  luminosity  of  the  compact  component  is 
L  =  (2. 3^4)  ’  10~4Lq  and  is  based  on  the  best  %2  value 
from  the  simultaneous  fit  of  Teff  and  R.  The  bounds  are 
the  minimum  and  maximum  luminosities  found  in  the  la 
contour  of  the  combined  x2-fit.  The  contour  levels  were 
derived  from  sorting  all  entries  in  the  likelihood  grid  and 
determining  those  values  of  the  likelihood,  where  the  cu¬ 
mulative  sum  up  to  these  values  correspond  to  certain 
confidence  levels  (e.g.,  la  corresponds  to  ~68%  confi¬ 
dence)  . 

In  Figure  5  we  show  the  spectral  energy  distribution 
of  the  best-fit  blackbody,  the  observed  fluxes  as  well  as 
a  set  of  representative  blackbody  curves  from  within  the 
ler  region. 

5.3.2.  Including  dust  extinction  effects 

In  case  the  object  is  embedded  in  the  mid-plane  of  the 
HD100546  circumstellar  disk,  the  dust  between  the  mid- 
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Fig.  5. —  Observed  flux  densities  (with  lcr  error  bars)  and  up¬ 
per  limit  (all  red  data  points)  over-plotted  on  blackbody  emission 
curves.  The  blue  line  corresponds  to  the  best  y2- fit  for  effective 
temperature  and  emitting  area.  The  blue  points  result  from  con¬ 
volving  the  blue  line  with  the  transmission  curves  of  the  Ks,  L' 
and  M'  filters.  The  black-dotted  lines  are  a  set  of  representative 
blackbodies  from  within  the  lcr  region  of  the  y2-fit. 

plane  and  the  disk  surface  layer  reduces  the  observed 
brightness  of  the  source  due  to  wavelength  dependent 
scattering  and  absorption  efficiencies  of  dust  grains.  We 
do  not  have  spatially  resolved  information  about  the  3D 
disk  structure  and  the  local  dust  grain  properties,  but  we 
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used  a  radiative  transfer  disk  model  for  HD100546  (Mul¬ 
ders  et  al.  2011;  Pineda  et  al.  2014)  to  estimate  the  opti¬ 
cal  depth  of  the  disk  at  the  location  of  the  protoplanet. 
The  values  were  2.52,  1.60  and  1.37  for  the  Ks ,  L'  and 
M'  filters,  respectively,  and  referred  to  the  optical  depth 
through  the  full  face-on  disk  at  ~50  au.  Using  these  es¬ 
timates,  the  observed  disk  inclination  and  assuming  that 
the  object  sits  in  the  disk  mid-plane  we  derived  extinction 
corrected  magnitudes  for  the  protoplanet  (see,  Table  2). 
With  the  extinction  corrected  values  for  the  magnitudes, 
we  re-ran  the  blackbody  fits  described  above  and  found 
an  effective  temperature  of  Teff  =  1242^^  K  and  a  ra¬ 
dius  of  the  emitting  area  of  R  =  7.3  ±3.2  Rjupiter-  These 
extinction-corrected  estimates  are  consistent  within  the 
error  bars  with  the  uncorrected  values  derived  in  the  pre¬ 
vious  section. 

We  emphasize  that  the  adopted  values  for  the  optical 
depth  are  entirely  based  on  a  circumstellar  disk  model 
and  with  the  data  in  hand  we  cannot  constrain  any  pos¬ 
sible  extinction  effects  empirically.  Furthermore,  simula¬ 
tions  suggest  that  the  formation  process  of  a  gas  giant 
planet  and  its  interaction  with  the  surrounding  circum¬ 
stellar  disk  is  a  three-dimensional  process  (e.g.,  Gressel 
et  al.  2013),  indicating  that  for  more  realistic  extinction 
estimates  other  effects  need  to  be  taken  into  account  as 
well.  However,  as  the  radiative  transfer  model  is  able  to 
reproduce  a  large  number  of  observational  constraints, 
we  consider  this  analysis  to  be  an  interesting  comparison 
with  the  default  results  based  on  the  observed  fluxes. 

5.4.  Comparison  with  evolutionary  and  atmospheric 
models  for  young  planets 

Given  the  available  data,  the  observed  morphology  and 
derived  parameters  are  best  explained  with  a  young,  po¬ 
tentially  still  forming,  gas  giant  planet  (cf.  Quanz  et  al. 
2013;  Currie  et  al.  2014).  Further  evidence  for  a  source 
orbiting  HD100546  around  50  au  comes  from  recent  ob¬ 
servations  with  the  Atacama  Large  Millimeter  Array 
(ALMA)  showing  that  the  distribution  of  mm-sized  dust 
grains  in  the  circumstellar  disk  mid-plane  hints  towards 
dynamical  interactions  with  a  young  protoplanet  (Pineda 
et  al.  2014;  Walsh  et  al.  2014). 

As  our  object  is  presumably  the  youngest  exoplanet 
discovered  so  far,  it  allows  for  a  direct  comparison  with 
model  predictions  for  the  earliest  stages  of  gas  giant 
planet  formation.  Models  with  high  values  for  the  ini¬ 
tial  entropy  (”hot-start“  models;  e.g.,  Marley  et  al.  2007; 
Baraffe  et  al.  2003)  predict  combinations  of  radius  and 
effective  temperature  that  agree  with  our  derived  param¬ 
eters  at  a  confidence  level  of  a  few  percent  or  even  less  in 
the  first  10  Myr  of  their  evolution  (see,  Figure  6).  The 
best  fit  is  found  for  an  object  with  ~5  Mjup;ter  at  an 
age  of  1  Myr.  For  older  ages,  the  best  fits  predict  higher 
masses.  In  general,  the  models  predict  smaller  radii  in 
the  relevant  temperature  range  (Baraffe  et  al.  2003;  Fort¬ 
ney  et  al.  2008).  As  we  will  discuss  in  section  6,  one  way 
to  explain  the  large  effective  radius  of  the  protoplanet  is 
to  assume  the  existence  of  a  spatially  unresolved  circum- 
planetary  disk  that  contributes  to  the  detected  flux.  All 
radius-temperature  fits  formally  improve  if  the  extinc¬ 
tion  effects  discussed  above  are  considered:  The  best-fit 
model  based  on  the  extinction-corrected  values  is  found 
for  a  1  Myr-old  ~10  MjUD;ter  object  with  a  confidence 
level  of  -10%  (cf.  Figure  6). 


X2  comparison  with  COND  models  - 
Simultaneous  fit  to  Teff  and  Radius 


Planet  mass  [MJupiter] 

Fig.  6. —  Quantitative  x2  comparison  of  results  from  our  black- 
body  fits  with  temperature-size  predictions  from  theoretical  models 
(Baraffe  et  al.  2003)  for  varying  masses  and  ages  of  young  gas  giant 
planets.  The  dash-dotted  line  shows  different  confidence  levels. 


X2  comparison  with  COND  models  - 
Simultaneous  fit  to  Ks,  L'  and  M'  fluxes 


Planet  mass  [MJupiler] 

Fig.  7. —  Results  from  a  x2-fit  of  the  predicted  Ks,  L ’  and  M' 
fluxes  of  young  gas  giant  planets  with  varying  masses  and  ages 
(Baraffe  et  al.  2003)  to  our  observed  fluxes. 

Going  one  step  further,  we  also  fitted  the  predicted  Ks, 
L'  and  M'  fluxes  from  the  COND  models  for  a  range  of 
masses  and  ages  of  1,  5  and  10  Myr.  It  is  important  to 
remember  that  the  predicted  fluxes  result  from  the  evo¬ 
lutionary  models  combined  with  additional  predictions 
from  the  atmospheric  models.  These  fits  yield  y2  values 
>25  (Figure  7),  which  is  significantly  worse  than  fitting 
radius  and  temperature  from  the  evolutionary  models 
alone.  However,  it  is  interesting  to  point  out  that  for 
a  given  age  the  best-fit  mass  range  is  basically  identical 
to  the  best-fit  mass  range  found  in  Figure  6.  In  Figure  8 
we  illustrate  the  main  differences  between  the  model  pre¬ 
dictions  and  the  observed  fluxes  for  those  regions  of  the 
mass-age  parameter  space  that  formally  provide  the  best 
fit  in  Figures  6  and  7.  It  shows  that,  irrespective  of  the 
assumed  age,  the  models  predict  Ks  magitudes  that  are 
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Fig.  8. —  Direct  comparison  of  Ks,  L '  and  M'  magnitudes  as  predicted  by  the  COND  models  (Baraffe  et  al.  2003)  with  the  observed 
magnitudes  of  the  protoplanet  for  ages  of  1,  5  and  10  Myr  (from  left  to  right).  The  grey-shaded  region  in  each  panel,  where  the  comparison 
is  shown,  is  defined  by  the  best-fit  mass  range  from  the  \2-fits  presented  in  Figures  6  and  7.  The  line  and  symbol  legend  shown  in  the  left 
panel  is  valid  for  all  panels. 


>3cr  discrepant  and  hence  should  have  led  to  a  detection. 
Furthermore,  in  all  cases  the  predicted  L'  magnitudes  are 
brighter  than  the  M'  magnitudes,  which  is  not  observed. 

Turning  to  models  with  low  values  for  the  initial  en¬ 
tropy  (” cold-star “  models),  typically  the  predicted  lumi¬ 
nosities  agree  with  our  derived  value  only  during  a  short 
phase  (<0.1  Myr)  at  the  beginning  and  at  the  end  of  the 
gas  runaway  accretion  (Marley  et  al.  2007).  This  would 
mean  we  have  caught  the  object  exactly  at  the  right  time, 
which  seems  unlikely.  More  recent  work  suggests  that  if 
the  solid  core  of  the  planet  consists  of  several  tens  of 
Earth  masses,  the  resulting  luminosity  might  be  compa¬ 
rable  to  what  we  found  here  even  a  few  million  years  af¬ 
ter  formation  (Mordasini  2013).  However,  the  predicted 
radii  of  the  objects  are  again  much  smaller  than  what  we 
derived  above. 

6.  DISCUSSION 

6.1.  The  protoplanet  HD100546  b 

One  of  the  most  interesting  questions  is  whether  we 
can  put  some  constraints  on  the  mass  of  the  protoplanet. 
While  all  ”  hot-start  “  models  predict  masses  of  at  least 
5  M jUpiter,  additional  observational  results  question  the 
presence  of  a  high-mass  protoplanet.  A  massive  planet 
orbiting  within  a  circumstellar  disk  opens  up  a  gap  with 
a  width  of  several  times  the  planets  Hill  radius  (>15  au  in 
this  case)  within  a  few  hundred  orbits  (Lin  &  Papaloizou 
1993).  However,  polarized  light  images  of  the  disk  show 
no  direct  evidence  for  a  disk  gap  (Avenhaus  et  al.  2014). 
As  the  orbital  time  of  the  object  is  only  ~250  years,  this 
would  imply  that  the  object  is  not  very  massive,  quite 
young  or  both.  Alternatively,  the  local  properties  of  the 
gas  disk  might  suppress  gap  formation  on  the  disk  surface 
(e.g.,  due  to  local  turbulent  viscosity)  or  the  combination 
of  disk  flaring  and  inclination  complicates  the  detection 
of  a  disk  gap  in  scattered  light.  Higher  spatial  resolu¬ 
tion  observations  in  the  future,  either  in  scattered  light 
or  with  ALMA,  can  help  us  to  search  for  clear  gap  signa¬ 
tures  at  the  object’s  location,  which  could  then  be  used 
to  put  some  constraints  on  the  object’s  mass. 

Concerning  the  discrepancies  between  the  observations 
and  the  model  predictions,  the  presence  of  a  circumplane- 
tary  disk,  as  expected  from  hydro-dynamical  simulations 
of  forming  gas  giant  planets,  helps  to  circumvent  some  of 


the  problems.  Such  a  disk  extends  out  to  roughly  40-50% 
of  the  planets  Hill  sphere  (e.g.,  Martin  &  Lubow  2011; 
Gressel  et  al.  2013)  and  would  add  an  additional  emis¬ 
sion  component  to  the  system.  Assuming  a  2  Mjup;ter 
protoplanet  its  circumplanetary  disk  would  be  ~1.4  au 
in  radius,  which  would  not  be  spatially  resolved  in  our 
images  and  would  hence  be  part  of  the  compact  emis¬ 
sion  component.  In  this  scenario,  our  derived  radius, 
but  of  course  also  the  derived  effective  temperature,  are 
then  the  superposition  of  the  emission  coming  from  the 
protoplanet  and  its  disk.  Such  a  two-component  model 
can  fit  the  data  for  several  reasonable  combinations  of 
sizes  and  effective  temperatures  and  data  points  at  ad¬ 
ditional  wavelengths  -  preferably  at  longer  wavelengths 
(Zhu  2015;  Eisner  2015)  -  are  needed  to  help  break  ex¬ 
isting  degeneracies.  In  principle  it  is  possible  to  artifi¬ 
cially  fix  the  effective  temperature,  mass  and  radius  of 
the  young  planet,  e.g.,  by  selecting  one  of  the  COND 
models,  and  then  use  the  observed  data  to  constrain  the 
properties  of  the  circumplanetary  disk.  However,  the 
COND  models  (just  like  any  other  evolutionary  model) 
are  highly  uncertain  and  unconstrained  by  empirical  data 
at  very  early  ages  and  furthermore,  given  the  results  in 
section  5.4.,  we  have  no  good  metric  to  decide,  which 
planet  model,  in  terms  of  age  and  mass,  we  should  pick. 

The  discovery  of  a  young  gas  giant  planet  still  embed¬ 
ded  in  the  circumstellar  disk  at  ~50  au  from  its  star 
suggests  that  these  objects  can  form  at  large  separa¬ 
tions.  This  is  in  particular  interesting  for  other,  slightly 
older,  exoplanet  systems,  where  massive  gas  giant  plan¬ 
ets  have  been  detected  at  comparable  orbital  separations 
(HR8799  be,  HD95086  b,  GJ504  b;  Marois  et  al.  2008; 
Rameau  et  al.  2013a, b;  Kuzuhara  et  al.  2013).  Also  some 
of  these  objects  may  have  formed  close  to  their  current 
location  and  additional  mechanism,  such  as  significant 
outward  migration  or  dynamical  scattering,  may  not  be 
required  to  explain  their  orbits.  However,  both  the  classi¬ 
cal  core  accretion  model  and  the  gravitational  instability 
model  cannot  easily  explain  the  data  presented  here.  In 
the  classical  core  accretion  model,  the  time  required  to 
build  up  a  rocky  core  of  several  Earth  masses  in  situ  at 
the  given  distance  from  the  star  exceeds  by  far  the  age  of 
the  system  (e.g.,  Kennedy  &  Kenyon  2008).  In  the  grav¬ 
itational  instability  model,  the  disk  has  to  be  massive 
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Fig.  9. —  Comparison  of  the  extended  emission  feature  southeast 
of  the  central  star  detected  in  our  ADI  data  with  the  location  of 
polarized  light  detected  in  polarimetric  differential  imaging  (PDI) 
data.  Specifically  we  compare  the  M'  ADI  image  (in  color),  where 
the  extended  emission  appears  to  be  more  pronounced,  with  the 
L'  PDI  image  from  Avenhaus  et  al.  (2014)  (contours). 

enough  to  locally  fragment.  The  remaining  mass  avail¬ 
able  in  the  HD100546  circumstellar  disk  (^10  Mjupjter, 
0.4%  of  the  stellar  mass;  Panic  et  al.  2010)  is  certainly 
not  sufficient  for  fragmentation  to  occur,  which  -  to  first 
order  -  sets  in  for  masses  of  ~10%  of  the  stellar  mass  (e.g., 
Lodato  2008).  So  even  if  the  disk  mass  was  a  factor  of 
a  few  higher  prior  to  the  formation  of  the  protoplanet, 
disk  fragmentation  seems  unlikely.  Recently  it  was  sug¬ 
gested  that  the  accretion  of  cm-sized  pebbles  that  are 
loosely  coupled  to  the  gas  in  the  circumstellar  disk  could 
significantly  speed  up  the  timescales  for  growing  a  rocky 
core  even  at  large  separation  from  the  central  star  (Lam- 
brechts  &  Johansen  2012).  Qualitatively,  such  a  model 
seems  to  be  able  to  explain  the  formation  of  a  gas  giant 
planet  at  ^50  au,  and  HD100546  might  be  an  ideal  lab¬ 
oratory  to  study  alternative  formation  processes  for  gas 
giant  planets. 

6.2.  The  extended  emission  features 

Finally,  concerning  the  extended  emission  component 
in  our  images,  this  is  likely  thermal  emission  from  warm 
material  in  the  surrounding  circumstellar  disk.  The  un¬ 
derlying  energy  source  is  unknown  at  the  moment,  and 
whether  or  not  local  compressional  heating  (e.g.,  Boley 
&  Durisen  2008),  or  similar  effects,  are  responsible  for 
the  observed  flux  requires  further  investigation.  As  de¬ 
scribed  in  section  4.2,  deriving  flux  estimates  from  our 
data  is  not  possible  without  substantial  modeling,  which 
is  beyond  the  scope  of  the  current  paper.  We  note  that 
in  our  final  images  another  region  of  extended  emission 
is  detected  to  the  southeast  of  the  star.  It  is  more  clearly 
seen  in  the  M'  images  but  is  also  present  in  L'  (see,  Fig¬ 
ure  1).  There  is  no  point  source  component  included 
in  this  emission  and  its  origin  is  also  unclear.  Part  of 
the  emission,  at  least  in  the  L'  band,  might  be  due  to 
scattered  light  from  the  disk  surface  as  in  this  direction 
from  the  star,  but  at  slightly  smaller  separations,  also 
the  L'  polarized  light  images  of  the  disk  surface  show  a 
flux  maximum  (Avenhaus  et  al.  2014).  We  show  a  di¬ 
rect  comparison  of  the  M'  image  presented  here  and  the 
polarized  light  image  in  Figure  9.  However,  given  the 
properties  of  typical  dust  grains,  the  scattering  efficiency 


and  hence  the  detected  flux  should  increase  with  shorter 
wavelength.  This  is  difficult  to  reconcile  with  our  detec¬ 
tion  in  the  M'  band  and  no  significant  emission  from  this 
structure  in  our  Ks  band  images.  The  feature  was  also 
reported  by  Currie  et  al.  (2014)  and  they  propose  that  it 
might  be  a  spiral  density  wave.  Indeed,  spiral  arm  fea¬ 
tures  have  been  reported  by  various  authors  and  at  var¬ 
ious  locations  in  the  HD100546  disk  (e.g.,  Grady  et  al. 
2001;  Ardila  et  al.  2007;  Boccaletti  et  al.  2013;  Avenhaus 
et  al.  2014).  More  data  is  required  to  determine  the  rela¬ 
tive  contribution  of  thermal  radiation  and  scattered  light 
to  the  observed  emission  and  to  understand  its  physical 
origin. 

7.  SUMMARY  &  CONCLUSIONS 

We  have  presented  the  first  multi-filter  study  of  the 
protoplanet  embedded  in  the  disk  of  the  Herbig  Ae/Be 
star  HD100546.  Our  key  results  can  be  summarized  as 
follows 

•  The  object  was  detected  in  L'  and  M'  and  con¬ 
sists  of  an  unresolved  point  source  component  and 
a  spatially  resolved  extended  emission  component. 
The  object  was  not  detected  in  Ks. 

•  The  contrast  of  the  point  source  component  relative 
to  the  host  star  is  9.4  ±  0.1  mag  and  9.2  ±  0.15 
mag  in  L'  and  M' ,  respectively.  The  3cr  limit  on 
the  minimum  contrast  in  Ks  is  9.6  mag.  These 
values  translate  into  apparent  magnitudes  of  L'  = 
13.92±0.10  mag,  M'  =  13.33±0.16  mag  and  Ks  > 
15.43  ±  0.06  mag. 

•  The  separation  between  the  point  source  compo¬ 
nent  and  the  host  star  is  (0.457  ±  0.014)"  and 
(0.472  ±  0.014)"  in  L'  and  M' ,  respectively.  The 
position  angle  is  (8.4  ±  1.4)°  and  (9.2  ±  1.4)°.  The 
average  de-pro jected  physical  separation  is  53  ±  2 
au. 

•  Combined  with  earlier  data  from  2011  we  demon¬ 
strated  that  the  object  is  co-moving  with  the  cen¬ 
tral  star,  and  also  the  L'  —  M'  color  and  apparent 
magnitudes  are  inconsistent  with  any  (sub-)stellar 
fore-  or  background  source.  Together  with  results 
from  other  studies  our  data  are  best  explained  with 
a  young  forming  planet  embedded  in  the  HD100546 
circumstellar  disk. 

•  Fitting  a  single  temperature  blackbody  to  the  ob¬ 
served  fluxes  of  the  point  source  component  yields 
an  effective  temperature  of  Teff  =  93213,02  K  and  a 
radius  for  the  emitting  area  of  R  =  6.9^2'g  R.Jupiter- 
The  best-fit  luminosity  is  L  =  (2. 3^0.4)  •  10 ~4L&. 
Teff  and  R  increase  when  possible  dust  extinction 
effects  caused  by  the  circumstellar  environment  are 
taken  into  account,  but  they  are  consistent  with  the 
values  above  at  the  ler  level. 

•  The  observed  L'  and  M'  magnitudes  are  inconsis¬ 
tent  (x2  >  25)  with  those  predicted  by  atmospheric 
models  for  young  (1-10  Myr)  gas  giant  planets. 
The  effective  temperature  and  radii  predicted  by 
evolutionary  models  for  young  gas  giant  planets 
agree  with  the  observations  at  the  <1%  level.  The 
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main  discrepancy  is  the  large  emitting  area  derived 
from  our  data. 

•  The  large  effective  emitting  area  of  the  object  can 
be  readily  explained  with  a  combination  of  a  young 
planet  and  a  surrounding  circumplanetary  disk.  In 
this  case,  the  derived  parameters  (Te f  f  and  R)  rep¬ 
resent  a  superposition  of  the  contributions  from 
both  components  (planet+disk)  as  the  circumplan¬ 
etary  disk  is  expected  to  be  unresolved  in  our  im¬ 
ages. 

Given  these  findings,  HD100546  is  a  unique  laboratory 
to  study  gas  giant  planet  formation  empirically.  Future 
ALMA  observations  with  comparable  resolution  as  the 
data  presented  here  will  confirm  the  existence  of  the  sus¬ 
pected  circumplanetary  disk  and  will  constrain  its  ex¬ 
tent  and  mass.  Such  observations  will  also  yield  spa¬ 
tially  resolved  information  about  the  physical  -  and  po¬ 
tentially  chemical  -  conditions  in  the  circumstellar  disk 
in  the  vicinity  of  the  forming  planet,  which  may  help  to 
further  constrain  the  processes  involved  in  the  object’s 
formation.  Finally,  new  high-contrast  imaging  observa¬ 
tions  with  VLT/SPHERE  or  Gemini/GPI  will  further 
push  the  detection  limits  at  Ks  or  even  shorter  wave¬ 
lengths  and  they  may  even  probe  directly  for  the  pre¬ 
dicted  planet  orbiting  within  the  disk  gap  at  ~  13-14  au 
(Brittain  et  al.  2014). 
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We  are  entering  an  era,  where  we  start  deriving 
empirical  constraints  on  the  formation  sites  and  for¬ 
mation  processes  of  gas  giant  planets,  and  together 
with  HD169142,  where  also  first  indications  for  multiple 
planet  formation  have  been  reported  (Reggiani  et  al. 
2014;  Osorio  et  al.  2014),  and  LkCal5  (Kraus  &  Ireland 
2012),  HD100546  will  be  one  of  the  prime  targets  for 
further  investigations. 
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